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Abstract 

An approach to the equation of state for the inner crust of neutron stars based on Skyrme- 
type forces is presented. Working within the Wigner-Seitz picture, the energy is calculated by the 
TETF (temperature-dependent extended Thomas-Fermi) method, with proton shell corrections 
added self-consistently by the Strutinsky-integral method. Using a Skyrme force that has been 
fitted to both neutron matter and to essentially all the nuclear mass data, we find strong proton 
shell effects: proton numbers Z = 50, 40 and 20 are the only values possible in the inner crust, 
assuming that nuclear equilibrium is maintained in the cooling neutron star right down to the 
ambient temperature. 

Convergence problems with the TETF expansion for the entropy, and our way of handling them, 
are discussed. Full TETF expressions for the specific heat of inhomogeneous nuclear matter are 
presented. Our treatment of the electron gas, including its specific heat, is essentially exact, and 
is described in detail. 



PACS numbers: 21.10.Dr, 21.60.Jz, 21.65+f, 26.60,+c 



I. INTRODUCTION 



We are concerned with the application of Skyrme-type effective nuclear interactions to 
the determination of the equation of state (EOS) of the inhomogeneous nuclear matter 
encountered at nuclear and subnuclear densities in core-collapse supernovas and in the inner 
crust of neutron stars. In our first paper on this topic |l( we adopted a Wigner-Seitz 
(WS) model of the inhomogeneous nuclear medium, and used the fourth-order semi-classical 
temperature-dependent extended Thomas-Fermi (TETF) method to calculate the kinetic 
energy and entropy. That paper dealt primarily with the conditions prevailing in core- 
collapse supernovas. The present paper relates rather to the inner crust of neutron stars, 
describing in particular some modifications to the earlier model, made necessary by two 
problems that emerge at the much lower temperatures T that are involved: i) the TETF 
expansion (in powers of h 2 ) for the entropy converges badly at low T; ii) proton shell effects 
are not negligible at low T. 

This last point is especially important if one is interested in the neutron-star crust as 
a possible alternative site for the synthesis of the so-called r-process elements 
The usual model of the r-process of nucleosynthesis is associated with the birth of a neutron 
star in a core-collapse supernova, during which "seed" nuclei are exposed to an intense 
flux of neutrons. Rapid ("r") capture of neutrons alternating with beta decay leads to the 
formation of a string of highly neutron-rich isotopes of a wide range of elements, which, 
once the source of neutrons is removed, will beta-decay back to the most neutron-rich stable 
isobar for the given mass number A (see Ref. 6[ for a recent review). The alternative 
picture, of interest here, is associated rather with the death of a neutron star, or at least 
with its partial disruption. Because of the very large densities, the matter in a neutron star 
is highly neutron rich, and the closer to the center the more neutron-rich it will be. But if for 
one reason or another matter is ejected from the neutron star it will rapidly decompress, and 
so will be able to undergo a chain of beta decays, the end-product of which will again be r- 
process nuclei. Ejection of matter from a neutron star is usually supposed to result from the 
merger of one neutron star with another, or with^a black hole 
been envisaged, e.g., volcanoes magnetars 

from the mass of the neutron star falling below the minimal critical value 
the precise ejection mechanism is of no concern to us in this paper. 
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but other scenarios have 



9|, quark stars [lOj and explosions resulting 



ll| . However, 



2 



It is convenient at this point to recall that at least three distinct regions can be recognized 
in a neutron star: a central, locally homogeneous, core, and two concentric shells charac- 
terized by different inhomogeneous phases [jjj]. The outermost of these shells, the "outer 
crust", consists of an electrically neutral lattice of nuclei and electrons. At the surface of 
the star only nuclei that are stable under natural terrestrial conditions are found (in fact, 
nuclear equilbrium, discussed below, implies that only 56 Fe will be found), but on moving 
towards the interior the increasing density leads to the appearance of nuclei that are more 
and more neutron rich, until at a mean local density p of around 2.4 xlO -4 nucleons.fm" 3 
(4.0 xlO 11 g.cm" 3 ) neutron drip sets in. This marks the transition to the "inner crust", 
which at least up to a mean density of p = 0.06 nucleons.fm -3 consists of neutron-proton 
clusters, or droplets, immersed in a neutron gas, with the neutralizing electron gas being 
essentially uniform (we neglect screening effects in this paper). It is equally well established 
that by the point where the mean density has risen to around p = 0.10 nucleons.fm -3 , i.e., 
about 2/3 of the density p of symmetric infinite nuclear matter (INM) at equilibrium, the 
droplet phase no longer exists and has been replaced by the homogeneous phase of the core, 
which consists primarily of neutrons, with a small admixture of proton-electron pairs, and 
possibly other particles, including free quarks, closer to the center. 

What happens in the transition region over the range 0.06 < p < 0.10 nucleons.fm" 3 , 
close to the inside edge of the inner crust, is far less clear. The question cannot be settled by 
observation at the present time, and theoretical predictions are sensitive to the details of the 
calculations, in particular to the choice of the effective interaction. For some interactions the 
transition from the droplet phase to the homogeneous phase is indirect and complex, with a 
whole sequence of different inhomogeneous phases being formed. At the interface with the 
homogeneous core these calculations find a "bubble" phase, this taking the form of bubbles 
of neutron gas in a denser liquid of neutrons and protons, the droplet phase having effectively 
been turned inside out. Furthermore, at slightly lower densities, between the bubble and 
droplet phases, several so-called "pasta" phases are predicted to put in an appearance, these 
being characterized by exotic, non-spherical shapes 121 ] . On the other hand, it has been 
shown that for other effective interactions the situation is much simpler, with no bubble or 
pasta phases being formed (at least at the assumed zero temperature of a stable neutron 
star): at a mean density of around p = 0.075 nucleons.fm -3 the droplet phase undergoes a 
transition directly to the homogeneous phase (see Ref. [131 ]. and references cited therein; also 
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Ref. [lj]). In the present paper paper we will avoid these ambiguities by limiting ourselves 
to values of p less than 0.06 nucleons.fm -3 , which means that we would not be able to deal 
with ejection mechanisms that reached even deeper into the star. 

Since neutron stars are formed at temperatures of the order of 10 MeV (10 11 K) and 



rapidly cool to around 0.1 MeV [12|, it is usually assumed that the final composition of 
the stable star corresponds to nuclear and beta equilibrium at a temperature of T = 0, the 
configuration of so-called "cold catalyzed matter" ; we shall later examine the validity of this 



assumption. Determining t 



re composition of the outer crust in this picture is straightforward 



see, for example, Ref. [15]): the equilibrating nucleus at each given density (or pressure) 



is found from the known nuclear masses, as given by experiment or, where mass data are 
unavailable, a mass model such as the FRDM [r| or HFB-14 [l^] (see also Refs. [3, 
for reviews). We shall therefore not consider the outer crust any further here. 

As for the composition of the inner crust of the stable neutron star, the relevant question 
at a given mean density p is to determine the total number of neutrons N, including those 
of the vapor, and protons Z per cluster. For this one needs the total Helmholtz free energy 
per nucleon / (including the electronic contribution) at the ambient temperature (usually 
taken to be zero), as a function of the density and the composition X = (Z, A = Z + N); 
one then minimizes / with respect to iV and Z at constant p. (Alternatively, to determine 
the composition at a given pressure P one minimizes the Gibbs free energy per nucleon g 
with respect to N and Z at constant P. It follows from the easily proven thermodynamical 
relation 



(9g\ = /df\ 



XI) 



\dX ) p,t \dX ) p,r 

that the two procedures are completely equivalent. We nevertheless find it more convenient 
to work with the Helmholtz free energy / at given values of p: see, for example, Section II 



of Ref. |20|.) 

The pressure in a layer of the crust of density p can then be found by numerical differen- 
tiation from the identity 



note particularly that / is a mean quantity, averaged over inhomogeneities, and not a local 
quantity. With the pressure P determined as a function of the mean density p, the values of 



P and p in any layer of the neutron star, along with the local composition, can be determined 



21 



22] 



through solution of the Tolman-Oppenheimer-Volkoff equation 

When, for one reason or another, decompression of crustal material begins, the temper- 
ature may start to rise. To follow the evolution of this process we shall require the EOS for 
non-zero values of T, and also the specific heat per nucleon at constant volume, c v , given in 
terms of the entropy per nucleon, s, by 



t(-) 

\dTJp,x 



;i.3) 



The entropy itself is given in terms of the Helmholtz free energy by 

_(d£\ 
S \dTJp,x ' 



;i.4) 



Thus all quantities of interest here can be derived from a calculation of / as a function of p 
and T. Note that s and c v , like /, are mean quantities. 

A popular EOS that has been extensively applied to supernova explosions is that of 



Lattimer and Swesty 



23l |. However, the applicability of this EOS to neutron-star crusts is 



limited by the fact that it is based on the so-called compressible liquid-drop model without 
any shell corrections, which at the prevailing low temperatures can be expected to be sig- 
nificant. Actually, both Refs. [2| and [,3f attempt to take account of shell effects, although 
in a rather rudimentary way, by making use of the algebraic bunching technique of Myers 
and Swiatecki (24]. 

In the present paper, as in Ref. I], we model the inhomogeneous nuclear medium by a 
single spherical WS cell, and attempt to incorporate shell effects into this framework micro- 
scopically and self-consistently, thereby permitting some measure of continuity of treatment 
across the interface between the inner and outer crusts. The most obvious way to do this is 
through the Hartree-Fock (HF) method, as has already been done, for example, by Bonche 



and Vautherin 



25( | at finite temperature and by Negele and Vautherin 26] at zero tempera- 



ture, using the WS approximation. However, we abandoned this approach for the following 
reason. While protons are strongly bound in the inner crust because of the large neutron 
excess, and thus show strong shell effects, for neutrons, by the very definition of the inner 
crust, there will be a continuous spectrum of unbound single-particle (s.p.) neutron states 
that are occupied. Thus any neutron added to the system must in reality go into this con- 
tinuum, whence it follows that we should not expect any neutron shell effects. Actually, 



this conclusion will hold only if the dripped neutrons form a uniform liquid, and in reality 
scattering of unbound neutrons on the mhomogeneities of the crust may give rise to so- 
called Casimir or band effects 
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29J|. whose exact evaluation requires the application 



of the band theory of solids (see Ref. 30j and references quoted therein). Nevertheless, 



these neutron shell effects are much smaller than the proton ones 



3l| and have therefore a 



negligible impact on the EOS and the equilibrium composition of the inner crust, although 



they are known to be significant for transport properties [32|. However, in practice any 
HF calculation in the WS approximation involves discretization, giving rise to shell effects 
for both protons and neutrons. But, as we have argued above, these neutron shell effects 



must be spurious, and in the HF calculation of Ref. 26j special steps had to be taken to 



smooth them (see also Refs. 



33(). We conclude that as far as neutrons are concerned 



the semiclassical extended Thomas-Fermi method is better adapted to a WS approach than 
is the HF method. 

The solution we adopt here to the problem of including the appropriate proton shell 
corrections without introducing spurious neutron shell corrections is to use the ETFSI (ex- 
ended Thomas- Fermi plus Strutinsky integral) high-speed approximation to the HF method 
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35. 36 



33, 



38]. We have already made an exploratory study of the applicability of this 



method to a WS picture of the EOS, and found that proton shell effects are indeed impor- 
tant [39j], but here, in addition to making much more extensive calculations of the EOS, we 
improve the T > results by taking account of possible shell effects in the entropy, which 
will manifest themselves in the free energy through the relation 



f = e-Ts 



;i.5) 



where e is the energy per nucleon. 

A further development of considerable significance is that the TETFSI method, as we 
shall refer to this temperature-dependent ETFSI method, is no longer limited to forces whose 
effective nucleon mass M* is equal to the real mass M. This permits us to use more realistic 
effective forces with smaller values of the effective mass. Thus in the present calculations 
the effective interaction that we use is the Skyrme force BSkl4, for which the effective 
mass in symmetric INM at the equilibrium density po (0.159 nucleons.fm -3 ) is 0.800M, 
which is to be compared with the value of 0.825M found in extended Brueckner-Hartree- 



Fock calculations that include three-nucleon forces 



401 ] . This is the force that underlies 
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the Hartree-Fock-Bogoliubov (HFB) mass model HFB-14 171 ]. a force that is eminently 
suitable for calculating the properties of neutron-star crustal matter, since on the one hand 
it has been fitted to the properties of neutron matter, as determined by calculations with 
realistic two- and three-nucleon forces [411 ]. and on the other hand it gives an excellent fit 
to essentially all the available mass data (<J rms = 0.729 MeV). Given that the neutron-star 
crust is both inhomogeneous and contains some protons, the high quality of the mass fit is 
especially relevant, since it means a) that inhomogeneities in nuclear matter (surface effects 
in droplet-model language) are well modeled, and b) that neutron-proton interactions are 
well represented. (However, no Skyrme force should be used for the highly supernuclear 
densities encountered deep within the core of a neutron star.) We stress that in this paper 



we neglect pairing, as in Refs. 



25 



26]. 



In Section II we discuss our parametrization of the WS cell. Section III describes our 
adaptation of the ETFSI method to the problem of the EOS of the neutron-star inner crust 
at non-zero temperatures, with particular attention to the convergence properties of the 
ETF expansion of the entropy. Our formalism is applied in Section IV to the properties of 
the inner crust of a neutron star (we do not examine in this paper the important question 
of the rapid decompression of neutron-star matter). The existence of strong proton-shell 
effects in the inner crust is discussed in this same section. Our conclusions are summarized 
in Section V. Some important material is to be found in the appendices, notably the TETF 
expansion of the specific heat (App. A) and a proof of the Strutinsky-integral theorem (App. 
C). 



II. THE WIGNER-SEITZ CELL 



Since we do not consider depths greater than that for which the mean density is p = 0.06 
nucleons. fm -3 only the droplet phase of nucleons has to be considered, and we shall assume 
that the WS cell associated with this phase is spherical, radius R c . This cell is entirely 
representative of the macroscopically sized volume element being considered, in the sense 
that all the nucleons of this volume element are imagined to be grouped into identical such 
cells, there being one cell for each droplet. The average neutron and proton densities over the 
locally representative cell must thus each be equal to the local values of the corresponding 
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macroscopic densities, p n and p p , given respectively by 

p n = p(l-Y e ) (2.1a) 

and 

P P = pY e , (2.1b) 

where Y e is the fraction of nucleons that are protons. The neutron and proton density 
distribution functions within the cell, p n (r) and p p (r), are then constrained by 



Lell d r R c JO 

where q denotes n or p, as the case may be, and the second equality holds in the case of 
spherical symmetry, assumed here. The total number of nucleons of each type in the cell is 

47T 

N q = -Rlp q (2.3) 

(N n = N,N p = Z). 

For the neutron and proton density distribution functions we adopt a modified version 
of the simple Fermi form that we used in Ref. [l|: limiting ourselves to the spherical case 
and writing 

Pq(r) = pBq + P0qfq(r) , (2.4) 

in which ps q is the usual constant background term, we now take 

f,(r) = !-2 -. . (2.5) 

C q —R c ^ 



Here the denominator of the radially varying term contains an extra factor ^ exp , ,, 
the presence of which guarantees that the density gradient now vanishes at the surface of the 
cell, thereby ensuring that the droplets merge smoothly with the homogeneous neutron vapor 
within the WS cell (in addition to being physically realistic, this condition is also required 
for the validity of the semi-classical part of the TETFSI method used here to calculate the 
nuclear kinetic energy and entropy). Also, again because of the extra factor, f q (r) itself 
vanishes on the surface of the WS cell, whence 

Pq(Rc) = PBq ■ (2-6) 



The factor 1/e ensures that the radially varying factor takes the value poq/2 at r = C q , as 
usual. Finally, we note that f q (r) varies monotonically over the cell. 

The two background constants ps q are not independent of the other parameters, but 
rather are fixed by 

3 

pBq = Pq- fflhPOq ' ( 2 ' 7 ) 

where 

h= 1^ f q {r)r 2 dr . (2.8) 
Jo 

For given p the WS cell is thus characterized by seven geometrical parameters, po q , C q , a q 
(q = n,p), and R c , in addition to the composition parameter Y e . In the most general case all 
eight of these parameters correspond to degrees of freedom, but either or both of the last two 
might be constrained to fixed values, according to the physical situation being described. In 
particular, for specified Z and N, the given value of p determines R c . 

But in all cases the complete set of parameters is subject to the additional constraint 
that the densities p q (r) must be positive at all points in the cell. We handle this problem as 
described in Section II of Ref. the modified radial distribution leading to the simplified 
condition 

P ' ! < Po q < % > (2-9) 



(3/,/ *g) - f q (r = 0) ^ 3I q 
given that f q {r) varies monotonically over the cell. The lower limit is seen from the mean- 
value theorem to be essentially negative, and since negative values of p 0q correspond to 
bubbles we shall be interested here only in the restricted range 

0<Po 9 <^ • (2.10) 

Although this is of no concern for the present paper, it is convenient to note here that 
once decompression and beta decay have begun, it will be necessary to constrain the back- 
ground parameters pB q to given values, since different (Z, A) configurations will be present 
simultaneously. Thus, with Eq. (12. 7p still holding, the two degrees of freedom corresponding 
to the parameters po q will be lost. Rather, they will be determined uniquely according to 

R 3 

POq = T^fipq - PBq) , (2.11) 



positive values still corresponding to droplets and negative to bubbles. The condition that 
the densities be everywhere positive then leads to constraints on the range of values that 
are possible for p Bq : 

°<p»< TZ 3V 4 /g(r - o )} • < 212 ' 

in which the upper limit is essentially larger than p q , but only the restricted range 

0<p Bg <p q (2.13) 

corresponds to droplets. 



III. THE TETFSI METHOD 



For a given set of the geometrical cell parameters we first write the total density of the 
Helmholtz free energy at a given point in the cell as 

T' = T nuc + F e + £ c + {p n M n + p p M p + n e m)c 2 , (3.1) 

where T nuc is the specifically nuclear free-energy density (discussed below), jF e is the density 
of the electron free energy (calculated essentially exactly), and £ c is the Coulomb energy. 
(Strictly speaking, in this equation p n and p p should be replaced by p n and p p , respectively, 
but the differences vanish on integrating over the cell.) Since p p = n e , electrical neutrality 
holding globally over the WS cell, we can now write 

(p n M n + p p M p + n e m)c 2 = p j (1 - Y e )M n + Y e (M p + m) }c 2 = pM n c 2 - pY e Q n<l3 , (3.2) 

where Q n ^ is the beta-decay energy of the neutron (0.782 MeV). But the term pM n c 2 makes 
a constant contribution to the free energy per nucleon and can thus be discarded. Thus in 
place of Eq. ( 13. lj) we write 

T = T nuc + + £c — pYeQn.fi ■ (3-3) 

The total Coulomb energy density, direct and exchange, is given in general by 

g c (r) = ^{ Pp {v) - n e } J - Q (pf + n^) , (3.4) 
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where for the exchange term we have used the usual Kohn-Sham variant [42[| of the Slater 
approximation. The last equation reduces in the case of spherical symmetry to 

p p {r')^—-r'^dr' + n e r -} (3.5) 

" sir 



The mean free energy per nucleon in the entire system is given by the corresponding 
quantity averaged over just one cell, 

f = jj H*)d 3 r , (3.6) 

where A = N + Z is the total number of nucleons in the cell. We shall likewise calculate the 
total entropy, the density of which at a given point can be written as 

a = cr nuc + a e = a n + cr p + a e , (3.7) 

whence for the mean total entropy per nucleon we have 

s = j I a(v)dh . (3.8) 



The Skyrme force BSkl4 
usual form 

1 



171 ] for which we calculate the densities T nuc and a nuc has the 



v 



'j 



t (l + x P (J )5(ry) + ti(l + x 1 P a )^{pi j 5{r ij ) + h.c.} 



+t 2 (l + x 2 Pa)^2Pij-S{r ij )pi j + ^t 3 (l + x 3 P a )p a 5{r ij ) 



1 

+ ^ (tr i + cr j ).p ij x% j )p ii . (3.9) 



The total nuclear energy density at any point can now be written as 

£- = E(^)+ v > ( 3 - 10 ) 

q 1 

where r q gives the kinetic-energy density of nucleons q as -^rT q (the first term of Eq. (3.10) 
is the kinetic-energy density multiplied by M q /M*), the effective mass M* is given by 

h 2 h 2 



2M* 2M q 



+ ^|t 1 (2 + x 1 )+t 2 (2 + x 2 )}p 



+ ^{t 2 (l + 2x 2 ) -t 1 (l + 2x 1 )}p q , (3.11) 
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and the static part of the potential energy by 

V = it |(2 + xo)p 2 -(l + 2x )^p^} 

g 

+ ^{3ti(2 + Xi) - t 2 (2 + x 2 )}(Vp) 2 

-^{3*i(l + 2x0 + t 2 (l + 2x 2 )} ^(Vp,) 2 

g 

+^s{(2 + x 3 )p 2 - (1 + 2x 3 ) Pl}p a 

g 

+ ^o$> q -V(p + p g ) , (3.12) 

g 

in which J q is the spin-current density, and we have now set the "quadratic current" term 
Vjj of Eq. (A4c) of Ref. [l| equal to throughout the calculation (our treatment of this 
term in Ref. |l( was inconsistent). For the nuclear free-energy density at any point we then 
have 

J r nuc = J2 IC + V > ( 3 - 13 ) 

g 

where 

h 2 

r=—-Ta-T(j Q . (3.14) 
g 2M* q 9 

The first stage of the full TETFSI method that we adopt in this paper for a given 
temperature T, mean density p and fixed values of N and Z consists in approximating the 
exact HF value of the nuclear free-energy density T nuc for the given Skyrme force by the 

43]: see, for example, our paper the appendix 



full fourth-order TETF method of BBD 



of which contains a convenient summary of the formalism, as we have used it here. We do 
not repeat this formalism here, although in Appendix [A] we present the TETF expansion 
for the specific heat, which appears not to have previously been published. Moreover, we 
remark here that the TETF expression for a q given by BBD [43J assumes that 

where it is to be noted that it is the local density p, and not the mean density p, that 
is held constant. We show in Appendix [B] that Eq. (13.151) follows from Eq. (I1.4p only if 



12 



equilibrium holds (this condition is made necessary by virtue of the temperature dependence 
of the density distribution itself). 

The essence of the TETF method is to express r q , J q and a q in terms of an assumed density 
distribution, which here we take to have the form given in Eqns. (12.41) and (12. 5p . The value 
that we obtain for the total free energy per nucleon /, as given by Eqns. (13.31) and (13.61) . 
is minimized with respect to the six geometrical parameters po q , C q , a q of the parametrized 
nucleon distribution. (In looking for full nuclear and beta equilibrium we minimize with 
respect to N and Z also.) The resulting nucleon distributions are denoted by p q , and the 
corresponding approximations to r q ,3 q ,a q and / by f q ,J q ,a q and /tetf respectively. The 
value of all these approximations to the exact HF values vary smoothly with respect to 
N and Z: it is a characteristic of the TETF method that shell corrections are lost. In 
the second stage we use the Strutinsky-integral method to correct J'tetf perturbatively for 
proton shell effects, as follows. 

At zero temperature, where /tetf = £etf, the corrected value of e takes the form 

e = e ETF + jE s p c , (3.16) 
where, according to the Strutinsky-integral theorem, 

dh(-^=,j p + p p U p + J p -W p ) . (3.17) 

2M* 

Here the integral goes over the volume of the WS cell, while the sum goes over all the 
occupied s.p. proton states, with the s.p. energies e^ p being the eigenvalues of the s.p. 
Schrodinger equation 

h 2 ~ ~ 1 ~ 

V + U p (r) - iWp(r) • V x a } <p Lp = e ijP (f) itP , (3.18) 

2M*(r) J 

and Hi the occupancy of s.p. state i (= or 1 for T = and no pairing). In these last 
two equations the effective mass M* is given by Eq. (13. lip with the smooth ETF densities 
p q replacing the exact HF densities p q . Likewise U p and W p are the central and spin- 
orbit proton fields, respectively, given by the usual HF expressions for these fields (see, for 



example, Eqns. (7) and (9) of Ref. [44]]), with p q replacing p q (note that here U p contains 
the Coulomb field). 

All three of these fields involve a folding of the Skyrme force over the nucleon distribution 
p q that emerges from the minimisation of cetf m the ETF part of the calculation. The fact 

13 



that the same Skyrme force underlies also the ETF part of the calculation implies a high 



degree of consistency between the two parts, which probably accounts 
found in comparisons with exact self-consistent HF calculations 



34 



or the close agreement 



351 ] (note that these tests 



were limited to the case of bound nuclei at T = and with M* = M q ). This theorem seems 
to have been stated for the first time in Ref. 



38 



45( | . A derivation, limited, however, to the 



case M* = M q , was sketched in Ref. 
effective mass, is presented in Appendix O 
As for T > 0, in Ref. 



. A more complete proof, applicable to arbitrary 



391 ] we replaced Eq. (13.161) by 



J^fTETF+jE^ , (3.19) 



where the proton shell correction is still given by Eq. (I3.17p . with 

n i = t— 77^ \~Fr\ ' ( 3 - 2 °) 

l + exp{(e i;P - Hp)/T} 

Hp being the chemical potential for protons. However, this does not take account of possible 
proton shell effects in the entropy, so here we will write rather 

fTETFSI = fTETF + jE S p C -T(s S /-- S T ETF ) , (3.21) 

where s*' p ' is the usual s.p. expression for the proton entropy, 



Kin ^ + (1-^)^(1-^)} , (3.22) 



in which the sum goes over all proton states. 

To determine the eigenvalues e^ p we expand the eigensolutions <pi yP to Eq. (13.181) in the 
basis defined by spherical Bessel functions ji(k n r) with the k n chosen to satisfy homogeneous 
boundary conditions (vanishing of the function or of its first derivative) on the surface of 
the WS cell. This generates a complete set of functions that are orthogonal over the cell, 
and we diagonalize the associated matrix. 

Interpolation Schemes. In addition to its much greater simplicity and rapidity, as com- 
pared to full-blown HF calculations, the (T)ETFSI method has the advantage of lending 
itself to interpolation. The point is that while the shell corrections are indeed fluctuating 
quantities, the method expresses these quantities in terms of quantities that themselves vary 
smoothly: the fluctuations arise entirely in the summation indicated on the right-hand side 
of Eq. (I3.17p . This feature was heavily exploited in the construction of the ETFSI-1 mass 
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table, the first mass table to be based on microscopic forces [38], and we anticipate that it 
could prove equally fruitful for the extensive tabulation of the EOS described in this paper. 

Specific heat. It is convenient to define a density of specific heat at constant volume, 
C v (r), at each point in the cell, according to 

Cv = jJ Cv{r ^ Y ' (3 ' 23) 

where c v is the specific heat per nucleon ( 11.31) . and the integration goes over the volume of 
the WS cell. It then follows, provided equilibrium holds, that 

(see Ap p endix IB"]) . This last result can then be used, starting from the TETF expansion 
of a [43j, to derive the corresponding expansion for C v up to order h 4 . This expansion is 
presented in Appendix [S] 

Convergence of the TETF expansions. 

Before applying the TETFSI method to the calculation of the neutron-star crust we make 
extensive tests of the convergence of the series expansion implicit in the TETF method. For 
this it is convenient to define the quantity 

4vr r ' 1 



t K q (r)r 2 dr , (3.25) 
A Jo 

which is just the "kinetic-thermal" part of the nuclear free energy per nucleon. 

We consider two densities, p = 0.06 and 3 xl0~ 4 nucleons.fm -3 , which lie close to the 
upper and lower limits, respectively, of the density range encountered in the inner crust. For 
the former we take a cell with Z = 30, A = 1123, and for the latter Z = 40, A = 150, these 
lying close to the equilibrium configuration (see Section HVT) . The results are shown in Tables 
HI and Hit respectively, for three temperatures, 0.1, 1.0, and 3.0 MeV. For each quantity, <fi n 
and <f) p , we show in successive lines the Thomas-Fermi (TF) approximation, the second- 
order correction in h, and the fourth-order correction; the sum of these three contributions 
is shown in the adjacent column. We also show in these tables the corresponding expansion 
terms for s n , s p , c„ jn , c V:P and the proton entropy as calculated by the s.p. expression (13.221) . 

We see that the expansions for all the neutron quantities, 4> n , (fi p , s n , and c v>ni converge 
well over the entire density and temperature range. However, the expansions for the proton 
entropy s p and the related specific heat c ViP diverge badly at low temperature, to the point of 
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occasionally giving negative entropies and specific heats. This proton-related problem is seen 
to be worse at low densities, which suggests that the fact that the neutron-related quantities 
converge so much better might be related to the fact that our neutron densities are always 
much higher than the proton densities. Indeed, we find that at an outer-crust density of p 
= 1 x 10~ 4 nucleons.fm -3 the TETF expansion for the neutron entropy of 208 Pb converges 
rather poorly. The existence of a significant background term, which gives no contribution to 
either the second- or fourth-order terms, appears to be crucial in this respect. It is fortunate 
that the lower density limit for good convergence of the TETF expansion for the neutron 
entropy lies below our domain of interest. In any case, it is quite clear that we cannot use 
the TETF expansions for the proton entropy or specific heat. 

For the proton entropy we therefore fall back instead on the s.p. expression (13.221) . which 
we require anyway for the proton shell corrections. On the other hand, we recall that we 
cannot use this s.p. expression for the neutron entropy because of continuum problems: it 
was precisely for this reason that we had to abandon the HF approach. As for the proton 
specific heat, we could in principle calculate numerically the temperature derivative of the 
proton entropy (13.221) and use Eq. (11.31) . but this is too time-consuming, if done accurately. 
We thus simply take the TF approximation for the proton contribution to the specific heat. 
Since the neutron contribution dominates the specific heat, the error thereby introduced will 
be relatively small; in any case we see that at higher temperatures the TF approximation 
to the entropy agrees approximately with the entropy calculated by Eq. ( 13.221) . 

It is reasonable to ask why the TETF entropy expansion diverges at low T. The second- 
order term contains a T~ l factor, and the fourth-order term a T~ 2 factor (see Eqns. (A. 25b- 
fl 

c) of Ref. [![), so that the corresponding numerators must likewise vanish at low T. This 
suggests that a massive cancellation within these numerators is leading to a significant loss 
of precision. However, we derived the low-T (strong degeneracy) limit of the TETF entropy 
expansion (see Appendix [D]) , and found essentially the same numerical results. We suggest 
that the reason for the observed breakdown in the TETF expansions of the entropy and 
specific heat lies with the fact that the validity of Eqns. ( 13. 15ft and ( 13.241) requires that the 
system be in equilibrium, which can never be exactly the case when the density profiles are 
parametrized, as in Eqns. (12. 4p and (12.51) . 
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IV. EQUILIBRIUM PROPERTIES OF INNER CRUST 



Pure TETF calculations. In calculating the equilibrium (nuclear and beta) composition 
of the inner crust, i.e., the number Z of protons and the number N of neutrons per WS 
cell at any given mean density p, we shall first neglect shell corrections. For this we simply 
minimize Jtetf, calculated at T = 0, with respect to N and Z. The results for Z, A = Z+N, 
and Y e = Z/A are shown in Figs. [U, [21 and[3l respectively. We stress that in this picture Z 
and A vary continuously, and are not restricted to integral values. 

These figures also show the values of the same quantities that we find for the force SLy4 



46( | . which has been widely applied to neutron stars. Fig. [2] for A shows that at a given p 



the WS cells tend to be bigger for force BSkl4, presumably because the value of the surface- 
energy coefficient a s f is slightly larger for the former (18.11 MeV as opposed to 17.6 MeV). 
However, since Y e runs slightly higher for SLy4 (see Fig. [3]) there are generally somewhat 
more protons for this latter force, as seen in Fig. [TJ (The higher values of Y e found for 
SLy4 can be traced to the fact that in homogeneous neutron matter the energy per neutron 
is higher for SLy4 over the entire density range considered here.) This sensitivity to the 
choice of force will have implications both for transport properties in the inner crust and for 
nucleosynthesis in decompressing neutron-star matter. 



Now the popular EOS of Ref. 



471 ] is also based on the SLy4 force, to the extent that the 



parameters of the underlying compressible liquid-drop model are calculated for this force. 
We compare the results of Ref. 47J for Z with our own in Fig. HI where it will be seen that 
the two models lead to appreciable differences even before taking shell effects into account. 

TETFSI calculations: proton shell effects. When shell corrections are included Z and A 
must be considered as integers, and thus change discontinuously: each (Z, A) pair defines 
a phase. Since p changes continuously transitions between one phase and another will take 
place over a finite range of p, in which the two phases coexist. We shall neglect this feature 
(see, however, the comments below), but our results will be valid over those considerable 
intervals of p for which only a single phase exists. 

The calculations proceed as outlined in Section IHIl Thus for given values of T and p, 
and for a given (Z, A) pair, fxETF has to be minimized with respect to the geometrical 
parameters of the cell, and then shell-corrected according to Eq. (13. 19[) (note that the 
(T)ETFSI method calculates shell corrections perturbatively). The equilibrium values of Z 
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and A, i.e., the values that minimize the TETFSI value of / for given values of T and p, are 
shown as functions of p for T = 0, 0.1, and 1.0 MeV in Tables UTTl HVl and|V] respectively. 

With the TETF values being shown in parentheses, we see in Tables IIHI and IIVI that at 
T = and 0.1 MeV there are strong proton shell effects, with Z = 50, 40, and 20 being 
successively favored as the density increases; in fact, these are the only values of Z that 
appear. The changes from Z = 50 at p = 0.005 fm -3 to Z = 40 at 0.01 fm -3 and from Z = 
40 at p = 0.04 fm~ 3 to Z = 20 at 0.05 fm~ 3 appear to be discontinuous. This situation can 
easily be understood from Figs. [5] and [6j where we show / at T = 0.1 MeV as a function of 
Z for p = 0.005 and 0.04 fm -3 , respectively, minimizing with respect to A for each value of 
Z: relatively strong minima occur for Z = 50, 40 and 20, and the system flips from one to 
the next as the density increases. (Somewhat surprisingly, although Z = 28 is consistently 
a strong local minimum it is never an absolute minimum.) It is noteworthy that the value 
of Z = 50 that we find here at the outside edge of the inner crust agrees with what we find 
for an outer-crust calculation using the HFB-14 mass model [ljj], i.e., a mass model based 
on the BSkl4 force used here. 

It will be seen from Figs. [5] and [6] that the minima are very close in energy, and numerical 
uncertainties in our computation often make it impossible to affirm with certainty which 
magic number prevails at a given density. Certainly, changing the force could be expected 
to lead to changes in the sequence of magic numbers. In fact, Negele and Vautherin 2(| 
report a quite difference sequence for the proton numbers 40 and 50 (and do not find Z = 
20 at any density). At T = 0.1 MeV the energy fluctuations associated with the Boltzmann 
factor represent an uncertainty on the total cell energy, and thus amount to 0.1/ A MeV per 
nucleon. Reference to Figs. [5] and [6] then shows that there will be negligible admixture of 
other values of Z with the magic values, although there could be significant admixture of 
A- values. 

Tables [TTT1 and HVI reveal no essential difference between T = and T = 0.1 MeV. However, 
with the energy differences between adjacent magic numbers being so small it could be that 
there are some intermediate densities for which the composition changes as T varies between 
and 0.1 MeV. This is one sense in which caution might have to be exercised in adopting 
the picture of "cold catalyzed matter". 

Table |V] shows that at T = 1 MeV shell effects have effectively been wiped out, even 
though there are still significant differences between the TETF and TETFSI values. This 



18 



raises the question of whether a nuclear equilibrium in the crust of a cooling neutron star 
can be maintained right down to T = 0.1 MeV. If "freeze-out" were to occur at, or slightly 
below, T = 1 MeV, i.e., if the complex rearrangement of nucleons necessary to maintain 
nuclear equilibrium could no longer take place over the lifetime of the neutron star, then 
the sharp shell effects that we have predicted here would not be observed. This is another, 
and probably more serious sense, in which the picture of "cold catalyzed matter" has to be 
carefully scrutinized. 

Pressure. We extract P directly from the computed values of frETFSi using Eq. (11.21) . 
numerically evaluating the derivative with a 3-point Savitzky-Golay filter (routine 'savgol' 



481]). Our results are shown in the penultimate columns of Tables [TTT | IIV } and El tests show 
that over this density range our results for P are reliable to within about 1 %. 

Another quantity of astrophysical interest is the temperature variation of the pressure 
|^ J ; using a well known Maxwell relation we have 

/ o 



This too requires a numerical differentiation, but only one: without using the Maxwell 
relation we would have to perform two numerical differentiations, with consequent loss of 
precision. Moreover, the derivative of s can be computed simultaneously with the one of / 
that gives us P, with negligible increase in computer time. The results are shown in the last 
columns of Tables Ifflt HV] and IVl 

Phase equilibrium. If we were to consider a quasi-continuum of values of p our calculations 
might show unphysical discontinuities in the pressure at the transition between the different 
(Z, A) phases. This is a result of our neglect of the possibility of a thermal equilibrium 
between the two phases in question, and in reality the pressure remains continuous. The 
transition pressure P is characterized by equality of the Gibbs free energy per nucleon in 
each phase, 

9i = 92 , (4.2) 

where 

9r = h + - ■ (4.3) 

Pi 

We stress that the equilibrium pressure satisfying this condition can be determined by calcu- 
lating fi as a function of p for each phase separately; in particular it is at no point necessary 
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to minimize the Gibbs function itself. Moreover, it can be shown that the condition (14. 2 j) 
follows from a minimization of the total Helmholtz function 49] . 



V. CONCLUSIONS 



We have developed here a high-speed approximation to the HF method for calculating 
the EOS of the neutron-star inner crust with Skyrme-type forces. Our method, which we 
refer to as the TETFSI method, models the inner crust in terms of the Wigner-Seitz cell, and 
consists essentially of a generalization to finite temperatures (and arbitrary effective mass) 



of the ETFSI method originally developed as a mass model [34], |35|, [36|, [37J, [38|]. An essential 
difference between our TETFSI method and a full-scale HF calculation of the EOS is that, 
whereas the latter method inevitably and automatically calculates both neutron and proton 
shell effects, here we calculate only the latter, since in reality shell effects are much weaker 
for neutrons than for protons, and will have negligible impact on the composition. In fact, if 
the HF method is used in a WS picture, as in the classical work of Negele and Vautherin 26|] 
it will lead, because of discretization, to spuriously large neutron shell effects As in Ref. 



26l | , we have neglected pairing in this paper, pending the determination of an effective pairing 



interaction appropriate to the conditions pertaining in neutron-star crusts. Nevertheless, it 
will be easy to inclu de p airing in the (T) ETFSI framework, as already done in the ETFSI 



mass models 



34 



35 



37|, 



38]. 



It was found that the TETF expansion of the entropy converges badly at low temperatures 
for densities typical of inner-crust protons (there was no problem for neutrons). We solved 
this difficulty by using the s.p. expression for the proton entropy. 

Our exploratory calculations of the EOS were performed with the Skyrme-type force 
BSkl4, a force that was fitted to essentially all the nuclear-mass data, forming thereby the 
basis of the HFB-14 mass model 12J. This force is particularly suitable for calculating the 
properties of neutron-star crustal matter, because it has been fitted to the properties of 
homogeneous neutron matter while at the same time the good fit to masses ensures that 
both inhomogeneities and the neutron-proton interaction are well represented. 

The calculated composition of the WS cells representing the clustering in the inner crust 
showed striking shell effects: for T = the proton number Z was limited to the magic 
values of 50, 40 and 20, the value decreasing with increasing density (at the interface with 
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the outer crust we found continuity with an outer-crust calculation based on the HFB-14 
mass model). Although essentially identical results are obtained for T = 0.1 MeV, all our 
calculated shell effects are wiped out at T = 1 MeV, which means that whether or not shell 
effects actually exist in the cold crust depends very much on the "freeze-out" temperature 
for nuclear equilibrium. On the other hand, we have shown that even without taking shell 
effects into account there are considerable differences between our predictions and those of 



471 ] is based. 



the compressible liquid-drop model on which the EOS of Ref. 

We intend to apply the method described here to a study of the synthesis of r-process 
nuclei in decompressing neutron-star crustal matter. To this end we present here, apparently 
for the first time, the TETF expressions for the specific heat of an inhomogeneous system of 
nucleons. In this same context of extensive computations over a wide range of temperature, 
density and composition, we point out that the (T)ETFSI method lends itself admirably to 
interpolation, without any loss of precision in the calculated shell effects, essentially because 
these arise in the sums of quantities that themselves vary smoothly 36]. 
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APPENDIX A: TETF EXPRESSIONS FOR THE SPECIFIC-HEAT DENSITY C v 

We expand the specific-heat density C v of Eq. (13.241) according to 



C V =C^+C^+C^ + --- , (A.l) 

where the first term on the right-hand side represents the Thomas-Fermi approximation, the 
second term the first-order correction in h 2 , and the last term the second-order correction 
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in h 2 . Then 



C5" = §{»f"+™-?*^} (A-*) 
Cfl = (24x 3 + 22x 2 + 5x - Q3xy - 33y + 45z) (A.2b) 

In Eq. (!A.2c|) we have 



G\ = (V2pg)2 (A.3a) 
Pg 

= v%(Vp,) 2 (A 3b) 

P. 

G5 = ^ (A.3c) 
P 3 , 



and 



%i\ f I- 1/2 (v q ) Y f 11 3 11 109 2 5 4 25 

— — = < ; — — > x xv x y H x H xz 

drj q J T \ I 1/2 (r) q ) I \192 60 y 320 y 64 64 

5 2 11 1 3 2 21 \ A . 

H x 2 H z y H y 2 w (A.4a) 

576 64 64 y 20 64 J v 1 

dX2\ f l-i/2(y g ) \ 2 ( 85 4 119 1 187 2 1 3 161 3 

t: — = < S — r- f x w H — z H v H x x 7/ 

\ h/ 2 (Vq) J \288 96 6 320 y 18 80 y 

609 2 175 159 2 H7 629 2 H 289 

+ xy xw H x z yz x y xy H xz 

320 y 64 64 64 y 480 y 60 y 192 

35 r 63 \ , . . 

H x 5 H v (A.4b) 

96 32 y v ; 

fdfcA f /-!/2fa) lV 77 11 4 33 2 693 391 3 2047 2 

— — = < —, — — } w H x H y u x y H xy 

\dr] q J T \ h/2(Vq) J V 256 192 256 y 256 240 y 1280 y 

719 4 161 97 o 63 2 1071 69 2 

— x y xw H x z x w H xv H x 2 

320 y 64 32 16 256 32 

207 315 11 2 11 161 - 483 

— yz H yw x y H xz H x H t> 

128 y 128 y 40 y 32 576 256 

21 , 4317 2 9 645 801 3 135 9 \ , k , 

+ 64* + 1280* V ~ U8 XyZ ~ Ym y + 128* J " (A ' 4c) 
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All quantities shown here are as defined in the Appendix of Ref. [l| and Appendix B of Ref. 
43) ]. except for u, which we define according to 

« = {h !f 2 • (A.5) 

Note that Eq. (1A.2c[) for the fourth-order term is valid only for an effective mass M* = M q 
(see the Appendix of Ref. 

As is the usual practice in expositions of the TETF method [43j, we have dropped terms 
here that vanish on integrating over configuration space, which in the present case corre- 
sponds to the WS cell. That is why we have required the density gradients to vanish at the 
surface of the cell (see Section HT1). 



APPENDIX B: PROOF OF EQNS. (13^ AND rf3^4l) 

The total free energy of the WS cell can be written as 

F = j F(p,Vp,T)d 3 r , (B.l) 
where the integration goes over the volume of the cell. For the entropy of the cell we have 

But, integrating by parts, we have 



,.~„ A®," 



\dVpJ 



^ V.(£Pl dh , (B.3) 



where we are making use of the vanishing of (dJ r /d'Vp) T on the surface of the cell; this 
follows from the fact that T must be at least quadratic in Vp, which must vanish on the 
surface of the cell for the TETF formalism to be valid. Then Eq. ( IB. 21) becomes 

--/©, A+ /(sm--(£)> • <B4) 
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Now if equilibrium holds at each temperature F must be a minimum with respect to varia- 
tions in p(r), which must therefore satisfy the Euler-Lagrange equation, 

where A is a Lagrange multiplier. Thus Eq. (IB .41) becomes 



But the second integral here is just the total number of nucleons in the cell (for simplicity 
we consider just one type of nucleon here), and since this is temperature independent Eq. 
0B.6p reduces to 

* = ' (B7) 

Eq. (13.151) follows at once. 

Likewise, for the specific heat we have from Eqns. (II .3p and (I3.2ip 

Then in exactly the same way as we have derived Eq. (IB. 41) . we find 

In the second term here we can write, using Eq. (I3.15p . 

(^) T " V - (av?). 



T 



But if we are at equilibrium Eq. (IB .51) will hold, and both sides of Eq. (IB. 101) will vanish, 
whence Eq. (IB.9P reduces to 



Ac v = T I d z r . (B.ll) 



Eq. (13.241) follows at once. 
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APPENDIX C: STRUTINSKY-INTEGRAL THEOREM 



To derive Eqns. (I3.16P and (13.171) we begin by noting that the exact HF energy of any 
finite nuclear system (nucleus or WS cell) for Skyrme forces and a Slater treatment of the 
Coulomb exchange energy can be written in the local form 

E HF = E HF [p,r,3] = y^{p(r),Vp(r),r(r),J(r)}rf 3 r , (C.l) 

where in terms of the exact HF s.p. functions 4>^ F (r) (not to be confused with the eigenso- 
lutions 4>i(r) of Eq. (I3.18P ) we have 

p(r) = 5>|0r(r)| 2 , (C.2a) 

i 

r(r) = ^n 4 |V0f F (r)| 2 (C.2b) 

i 

and 

J(r) = -^n J 0f F *(r)Vx < T0r(r) ; (C.2c) 

i 

for simplicity we do not distinguish here between the two charge states. For the exact HF 
quantities p(r),r(r) and J(r) let us now write 

p = p + Sp , (C.3a) 



t = t + 5t (C.3b) 

and 

3 = 3 + 63 , (C.3c) 

where p, f and J represent the smooth quantities emerging from the (T)ETF calculation. 
Then to first order in Sp, St and S3 we have 

E H f = J £{p(r),Vp(r),.-. , f(r), J(r)}rf 3 r 

+ /{(SL^(iL-(iL-"} A - 
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where for the functional derivatives appearing here we have 

'5£ 



8p 



= U(r) , (C.5a) 

p,f,3 



and 



-2 



Pi =^t- ( c.6b) 

*rJ?,Tj 2M*(r) 



g) = W W , (C,c) 



in which [/, ?i /(2M*) and W are the smoothed ETF fields appearing in Eq. (13. 18ft . Thus 
Eq. ( 1C.4I) becomes 



/:;/// = E E t\- I < U(r)Sp + -JL — St + W(r) • 53 I r/ :! r . (C.6) 

2M*(r) J 

Comparing with Eq. (13.161) . we see that we now have to identify the integral in this last 
equation with the shell correction E sc of Eq. ( 13.17ft . 

We next replace the exact HF s.p. functions <fif F (r) in Eqns. (IC.2aj) . (IC.2bl) and (IC.2cj) 



by the eigensolutions 0,(r) to Eq. (I3.18p . and define thereby the quantities p',r' and J', 
respectively; these quantities will certainly be fluctuating. We can then write 

5p = ( p '-p) + (p-p>) , (C.7) 

and likewise for 5t and 53. We recall now that p is the initial approximation (ETF) to the 
exact HF density p, while p' represents our attempt to improve on this approximation. Thus 
the first term on the right-hand side of Eq. (IC.7I) can be regarded as the first-order estimate 
of the correction to p, while the second term represents the residual error. Accepting, then, 
that our method is essentially one of first-order perturbation theory, we simply drop the 
second term of Eq. ( 1C.7I) . Treating 5t and 53 in the same way, Eq. (1C.6I) reduces to 

E HF = E ETF + / J U(r)(p' -p) + -J^—(t' -t)+ W(r) • (J' - J) I dh . (C.8) 
J \ 2M*(r) J 

But it follows from Eq. (13TT81 that 

/ \u(r)p> + J^r' + W(r). 3] dh = J2n^ • (C.9) 
J ( 2M*(r) J V 
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Then Eq. (1C.8j) reduces to 

E HF = E ETF + V - ( | U{r)p + -J^-f + W(r) • J) I d 3 r . (CIO) 
i J l 2M*(r) J 

This completes the proof that to first order the shell correction is given by Eq. (13.171) . 

APPENDIX D: STRONG-DEGENERACY LIMIT OF THE TETF EXPANSION 
FOR ENTROPY 

In the limit of low temperature and high density the TETF expansion of the entropy 
density a q is as follows. For the Thomas-Fermi approximation we have 

The first-order correction in h 2 to this is 

( 2 ) vr 2 M g 1 1 \ {V Pq f 9 fV/,\ 2 3 



* " 27 ft 2 (37r 2 p 9 )4/3 /g | pq ^ \ fq J fq 

+ 0(T 3 ) , (D.lb) 

while the second-order correction in h 2 is 

tt 2 M n 1 / „ 413 



^ = —r—^-, (l7G? - ^Gg + fo^j T + 0(T 3 ) . (D.lc) 



q 1620 h 2 (3ir 2 p q ) 2 



In these equations we have defined f q = M q /M*, while the G\ are defined in Eqns. (IA.3al) . 
flA73b|) and flA3cD . 
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TABLE I: Convergence of the TETF expansion at p = 0.06 nucleons.fm 3 for cell with Z = 30, A = 
1123. Temperature T in MeV. 



T 


0.1 


1.0 


3.0 




0.1943E+02 
-0.1027E-01 19.42 
0.3566E-03 


0.1934E+02 
-0.9452E-02 19.33 
0.2959E-03 


0.1845E+02 
-0.3529E-02 18.45 
0.8834E-04 


4>p 


0.1818 
0.3598E-02 0.1877 
0.2299E-02 


0.1659E+00 
0.5124E-02 0.1724 
0.1404E-02 


-0.2328E-01 
0.2516E-02 -0.02054 
0.2167E-03 


Sn 


0.1477E-01 

-0.1387E-05 0.01477 
-0.3124E-07 


0.1476E+00 
-0.1318E-04 0.1476 
-0.2757E-06 


0.4401E+00 
-0.1879E-04 0.4401 
-0.3389E-06 


Sp 


0.2110E-02 

-0.1537E-02 -1.005E-03 
-0.1578E-02 


0.1695E-01 
-0.2133E-02 0.01516 
0.3437E-03 


0.6378E-01 
-0.4283E-03 0.06337 
0.2441E-04 


s.p. 
Sp 


0.0036233 


0.01646 


0.06446 


Cv,n 


0.1477E-01 

-0.1388E-05 0.01477 
-0.3125E-07 


0.1473E+00 
-0.1352E-04 0.1473 
-0.2880E-06 


0.4316E+00 
-0.2503E-04 0.4316 
-0.5186E-06 


Cv,p 


0.1895E-02 
-0.9657E-03 0.0075 
0.6571E-02 


0.1372E-01 
0.4703E-03 0.01379 
-0.4022E-03 


0.3366E-01 
0.4400E-03 0.03410 
-0.2655E-05 
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TABLE II: Convergence of the TETF expansion at p = 3 x 10 4 nucleons.fm 3 for cell with 
Z = 40, A = 150. Temperature T in MeV. 



T 


0.1 


1.0 


3.0 




0.1472E+02 
-0.2021E+00 14.57 
0.4855E-01 


0.1164E+02 
-0.1824E+00 11.52 
0.5564E-01 


-0.2114E+01 
-0.5707E-01 -2.135 
0.3647E-01 


(f>p 


0.4465E+01 
-0.3411E-01 4.496 
0.6509E-01 


0.4699E+01 
-0.3642E-01 4.727 
0.6379E-01 


0.6777E+00 
-0.1831E-01 0.6960 
0.3655E-01 


s n 


0.2082E+00 
-0.8479E-02 0.1936 
-0.6098E-02 


0.1026E+01 
-0.2113E-01 1.006 
0.1083E-02 


0.2923E+01 
-0.1382E-01 2.909 
0.1035E-03 


Sp 


0.7193E-02 

-0.1315E-01 0.01024 
0.1620E-01 


0.6150E-01 
-0.2234E-01 0.0405 
0.1303E-02 


0.8998E+00 
-0.1353E-01 0.8864 
0.1074E-03 


s.p. 
Sp 


0.1199E-03 


0.06322 


0.7088 


Cv,n 


0.1385E+00 
-0.6991E-02 0.1410 
0.9451E-02 


0.3958E+00 
-0.4233E-02 0.3905 
-0.1176E-02 


0.8998E+00 
0.1726E-02 0.9015 
-0.1625E-03 


Cy,p 


0.6965E-02 

-0.4111E-02 -0.00867 
-0.1152E-01 


0.5660E-01 
-0.2989E-02 0.05094 
-0.2666E-02 


0.2593E+00 
0.2649E-02 0.2619 
-0.3149E-04 
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TABLE III: TETFSI results for number of protons Z and total number of nucleons A in WS cell 
for nuclear and beta equilibrium at T = as a function of p for force BSK14. TETF results in 
parentheses. Last two columns show TETFSI values of pressure P and ( ) _ 



p (fm 3 ) 


Z A 


P (MeV.fm- 3 ) 


(f)-(f-- 3 ) 


0.0003 


50 (38) 200 (146) 


0.000940 





0.001 


50 (39) 460 (385) 


0.00179 





0.005 


50 (39) 1140 (831) 


0.00813 





0.01 


40 (38) 1215 (1115) 


0.0185 





0.02 


40 (35) 1485 (1302) 


0.0448 





0.03 


40 (33) 1590 (1303) 


0.0784 





0.04 


40 (31) 1610 (1261) 


0.121 





0.05 


20 (30) 800 (1171) 


0.175 





0.06 


20 (29) 780(1105) 


0.243 
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TABLE IV: As for Table HI] but with T = 0.1 MeV. 



p (fm 3 ) 




Z 


A 


P (MeV.fm- 3 ) 




0.0003 


50 


(38) 


200 (147) 


0.000946 


0.0000423 


0.001 


50 


(39) 


460 (341) 


0.00179 


0.000140 


0.005 


50 


(38) 


1130 (842) 


0.00816 


0.000270 


0.01 


40 


(38) 


1210 (1107) 


0.0185 


0.000346 


0.02 


40 


(35) 


1480 (1294) 


0.0448 


0.000444 


0.03 


40 


(33) 


1595 (1303) 


0.0784 


0.000511 


0.04 


40 


(31) 


1610 (1242) 


0.121 


0.000568 


0.05 


20 


(30) 


800 (1190) 


0.175 


0.000617 


0.06 


20 


(29) 


765 (1116) 


0.243 


0.000662 



TABLE V: As for Table HI] but with T = 1.0 MeV. 



p (fm 3 ) 




Z 


A 


P (MeV.fm- 3 ) 


(§)_ (f-- 3 ) 


0.0003 


46 


(37) 


310 (234) 


0.000633 


0.000182 


0.001 


46 


(38) 


520 (450) 


0.00192 


0.000631 


0.005 


44 


(39) 


1020 (858) 


0.00936 


0.00233 


0.01 


42 


(37) 


1280 (1120) 


0.0202 


0.00329 


0.02 


40 


(36) 


1480 (1307) 


0.04701 


0.00434 


0.03 


38 


(33) 


1505 (1301) 


0.0810 


0.00502 


0.04 


36 


(31) 


1450 (1232) 


0.124 


0.00553 


0.05 


34 


(30) 


1340 (1165) 


0.179 


0.00568 


0.06 


26 


(29) 


985 (1082) 


0.246 


0.00445 
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FIG. 1: Number of protons Z in WS cell given by ETF method for nuclear and beta equilibrium 
at T = as a function of density (fm -3 ) for forces BSK14 and SLy4. 
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FIG. 2: Total number of nucleons A = Z + N in WS cell given by ETF method for nuclear and 
beta equilibrium at T = as a function of density (fm -3 ) for forces BSK14 and SLy4. 
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FIG. 3: Fraction Y e = Z/A of nucleons that are protons given by ETF method for nuclear and 
beta equilibrium at T = as a function of density (fm~ 3 ) for forces BSK14 and SLy4. 
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FIG. 4: Comparison of ETF and CLD (compressible liquid drop) calculations of equilibrium value 
of Z at T = as a function of density (fm~ 3 ) with SLy4 force. 
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FIG. 5: Variation of /t et fsi with Z (always for optimal value of A) at p = 0.005 fm 3 and T 
0.1 MeV. 
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FIG. 6: Variation of Jtetfsi with Z (always for optimal value of A) at p = 0.04 fm 3 . and T 
0.1 MeV. 
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